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Abstract 

We consider some (anisotropic and piecewise constant) convection-diffusion-rea- 
ction problems in domains of M 2 , approximated by a discontinuous Galerkin method 
with polynomials of any degree. We propose two a posteriori error estimators based 
on gradient recovery by averaging. It is shown that these estimators give rise to an 
upper bound where the constant is explicitly known up to some additional terms that 
guarantees reliability. The lower bound is also established, one being robust when the 
convection term (or the reaction term) becomes dominant. Moreover, the estimator is 
asymptotically exact when the recovered gradient is superconvergent. The reliability 
and efficiency of the proposed estimators are confirmed by some numerical tests. 

Key Words convection-diffusion-reaction problems, a posteriori estimator, discontinuous 
Galerkin finite elements. 
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1 Introduction 



The finite element method is the more popular one that is commonly used in the numerical 
realization of different problems appearing in engineering applications, like the Laplace 
equation, the Lame system, the Stokes system, the Maxwell system, etc.... (see [?] |8| 126]). 
More recently discontinuous Galerkin finite element methods become very attractive since 
they present some advantages. For example, they allow to perform "p refinement", by 
locally increasing the polynomial degree of the approximation if needed. They can moreover 
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use non-conform meshes allowing hanging-nodes, making the mesh generation easier for 
concrete industrial applications. We refer to (3J [10], and the references cited there, for a 
good overview on this topic. Adaptive techniques based on a posteriori error estimators 
have become indispensable tools for such methods. For continuous Galerkin finite element 
methods, there now exists a vast amount of literature on a posteriori error estimation for 
problems in mechanics or electromagnetism and obtaining locally defined a posteriori error 
estimates, see for instance the monographs [21 21 [271 ESj • On the other hand a similar theory 
for discontinuous methods is less developed, let us quote (51 [121 [131 I2Q1 I2H 1221 [301 EI]- 

Usually upper and lower bounds are proved in order to guarantee the reliability and 
the efficiency of the proposed estimator. Most of the existing approaches involve constants 
depending on the shape regularity of the elements and/or of the jumps in the coefficients; 
but these dependences are often not given. Only a small number of approaches gives rise 
to estimates with explicit constants, let us quote [21 El [HI [231 ESI (23 129] for continuous 
methods. For discontinuous methods, we may cite the recent papers [H [HI EH E2 EE [21] • 

Our goal is here to consider convection-diffusion-reaction problems with discontinuous 
diffusion coefficients in two-dimensional domains with Dirichlet boundary conditions ap- 
proximated by a discontinuous Galerkin method with polynomials of any degree. Inspired 
from the paper [18] , which treats the case of continuous diffusion coefficients approximated 
by a continuous Galerkin method, we further derive some a posteriori estimators with an 
explicit constant in the upper bound (1 or a similar constant) up to some additional terms 
that are usually superconvergent and some oscillating terms. The approach, called gradient 
recovery by averaging [18] is based on the construction of a Zienkiewicz/Zhu estimator, 
namely the difference in an appropriate norm of aVhUh ~ Gut, where VhUh is the broken 
gradient of Uh and Guh is a H(div )-conforming approximation of this variable. Here special 
attention has to be paid due to the assumption that a may be discontinuous. Moreover the 
non conforming part of the error is managed using a comparison principle from [TH] and 
a standard Oswald interpolation operator [H [22]- Furthermore using standard inverse in- 
equalities, we show that our estimator is locally efficient. Two interests of this approach are 
first the simplicity of the construction of Guh, and secondly its super convergence property 
(validated by numerical tests). Let us finally notice that this paper extends our previous 
one [13] in many directions: first we treat convection-diffusion-reaction problems instead of 
purely diffusion ones. Second we track the dependence of the constant in the lower bounds, 
in particular we show that the natural extension of the estimator from [13] yields a lower 
bound with a constant that is not robust when the convection and/or the reaction term 
becomes dominant. Consequently we introduce another estimator (adapted from [36J) that 
is robust, keeping nevertheless an upper bound with an explicit constant. 

The schedule of the paper is as follows: We recall in section 2 the convection-diffusion- 
reaction problem, its numerical approximation and recall the comparison principle from 
|16j . Section 3 is devoted to the introduction of the first estimator based on gradient 
averaging and the proofs of the upper and lower bounds. The upper bound directly follows 
from the construction of the estimator and some results from [T8] , while the lower bound 
requires the use of some inverse inequalities and a special construction of Gu^. Since the 
lower bound is not robust as the convection and/or the reaction term becomes dominant we 
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propose an alternative estimator in section 4 and show its robustness. Finally in section 
[5] some numerical tests are presented that confirm the reliability and efficiency of our 
estimators and the superconvergence of Guh to dWu. 

Let us finish this introduction with some notation used in the remainder of the paper: 
On D, the L 2 (D)-norm will be denoted by || ■ ||/> In the case D = Q, we will drop 
the index Q. The usual norm and semi-norm of H S (D) (s > 0) are denoted by || ■ \\ St o 
and | • \ St D, respectively. Finally, the notation a < b and a ~ b means the existence of 
positive constants C\ and C*2, which are independent of the mesh size, the coefficients of 
the operator and of the quantities a and b under consideration such that a < and 
Cib < a < respectively. In other words, the constants may depend on the aspect ratio 
of the mesh, as well as the polynomial degree /, but they do not depend on the coefficients 
of the operator a, (3 and /1 (see below). 



2 The boundary value problem and its discretization 

Let Q be a bounded domain of M 2 with a Lipschitz boundary T that we suppose to be polyg- 
onal. We consider the following convection-diffusion-reaction problem with homogeneous 
Dirichlet boundary conditions: 

J — div (a Vw) + (3 ■ V« + jiu — f in Q, , . 

\ u = on T. ^ > 

In the sequel, we suppose that a is a symmetric positive definite matrix which is piece- 
wise constant, namely we assume that there exists a partition V of Q into a finite set 
of Lipschitz polygonal domains Qi, ■ ■ ■ ,Qj such that, on each Qj, a = aj where aj is a 
symmetric positive definite matrix. Furthermore we assume that (3 G if (div , Q) nL°°(f2) 2 , 
/i G L°°(Q) and are such that /i — |div (3 > 0. If [i — |div (3 = on cu C f2, then we assume 
that /i = div (3 = on uj. 

The variational formulation of ([T|) involves the bilinear form 

B(u,v)= / (aVu- Vv + f3- Vuv + fmv),Vu,v G Hq(SI), 
Jn 

and the corresponding energy norm 

IIMH 2 = B{v,v) = I (aVv ■ Vv + (// - ldiv(3)\v\ 2 ), Wv G H%{Q). 

Given / G L 2 (Q), the weak formulation consists in finding u G -^o(^) such that 

B(u,v) = (f,v), WveH^Q), (2) 

where (/, v ) means the L 2 inner product in Q, i.e., (/, v) = j n fv. Invoking the positiveness 
of a and the hypothesis (fi — |div 0) > 0, the bilinear form B is coercive on Hq(Q) with 
respect to the norm 1 1 1 • 1 1 1 and this coerciveness guarantees that problem ([2]) has a unique 
solution by the Lax-Milgram lemma. 
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2.1 Discontinuous Galerkin approximated problem 

Following [31 [161 [22], we consider the following discontinuous Galerkin approximation of 
our continuous problem: We consider a triangulation T made of triangles T whose edges 
are denoted by e. We assume that this triangulation is regular, i.e., for any element T, 

the ratio — is bounded by a constant o > independent of T and of the mesh size 
Pt 

h = maxx^T hr, where hr is the diameter of T and px the diameter of its largest inscribed 
ball. We further assume that T is conforming with the partition V of Q, i.e., the matrix a 
being constant on each T £ T, we then denote by ax the value of a restricted to an element 
T. With each edge e of an element T, we associate its length h e and a unit normal vector 
n e , while stands for the outer unit normal vector along dT. E (resp. Af) represents 
the set of edges (resp. vertices) of the triangulation. In the sequel, we need to distinguish 
between edges (resp. vertices) included into Q or into T, in other words, we set 

E int = {e £ E : e C n},M int = {x £ Af : x £ Q}, 
E D = {e £ E : e C T},Af D = {x e Af : x £ T}. 

Problem ()2]) is approximated by the (discontinuous) finite element space: 

x h = {v h £ L 2 (n)\v h]T £ p,(t),t £ r} , 

where £ is a fixed positive integer. Later on we also need the continuous counterpart of 
Xh, namely we introduce 

S h = {v h £ C(H)K| T £ F e {T),T £ T} , 

as well as 

S hA = {v h £ C(fi)|v ft | T £ Pi(T),T £ T} . 

We further need 

X M = {v h £ L 2 (fi)|^| T £ P!(T),T £ T} . 

For our further analysis we need to define some jumps and means through any edge 
e £ E of the triangulation. For e £ E int , we denote by T + and T _ the two elements of T 
containing e. Let g £ X^, we denote by q ± , the traces of q taken from T ± , respectively. 
Then we define the mean of q on e by 

{{<?}L = w + (e)g + + u-(e)q-, 

where the nonnegative weights w ± (e) have to satisfy u + (e)+u~(e) = 1. If w + (e) = w~(e) = 
1/2, we drop the index uj. The jump of q on e is now defined as follows: 

[q] = q + n T + + q~n T -. 

Remark that the jump \q\ of q is vector-valued. 
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For v G [Xh] , we denote similarly 

i»L = w + (e)i; + + w"(e)v". 

For a boundary edge e, i. e. e G £d, there exists a unique element T + G T such that 
e C dT + . Therefore the mean and jump of q are defined by {{q}} = q + and \q\ = q + riT+- 
For q G Xh, we define its broken gradient V^g in Q by : 

(V&g)| T = V^VTgT. 

For io G X^ + i^o(^) aric ^ f° r w C O, we set 

/ (aVw ■ Vw + (fi - ^div (3)\w\ 2 ),VT £ T, 
Jt 2 

E iihiIt, \\h\\ 2 = J2 



w 



|2 

It 



FIIL = 2^ IIFNIt> IIFIII = 2^ HFI 

|2 



1 2 



H\DG,h = \\\W\\\ + [J2 h ^\ 



Note that || • \\DG,h is a norm on Xh. 

With these notations, we define the bilinear form Bh(., .) as follows: 

u ft ) = (aV h u h , V h v h ) + ((/i - div (3)u h , v h ) - (u h , f3 ■ V h v h ) (3) 
- E / i d {{<^hV h }} u ■ K] + {{aV fc u fc }} w • KD 

+ E / (^M ' M + {{uh}}? ■ M), Vu ft , g x h , 

where 9 is a fixed real parameter and the positive parameters 7 e are chosen appropriately, 
see below. 

The discontinuous Galerkin approximation of problem (J2J) reads now: Find Uh G Xh, 
such that 

B h {u h , v h ) = (f, v h ), Vv h G X h . (4) 

In taking the interior weights w ± (e) equal to 1/2 and setting = 0, 9 = — 1 or 
9 = 1 leads to the incomplete, nonsymmetric or symmetric interior-penalty discontinuous 
Galerkin methods. The stabilization parameters 7 e are chosen in the form (see e.g. [ISJdT]) 

7e = a e-T— + 7/3,e, (5) 
h e 

where a e is a positive parameter, 7 a e is a positive parameter that depends on a and 7^ 
is a positive parameter that depends on (3 and is zero is (3 = (the choice 7^ = 
corresponding to an upwinding scheme). Whenever 9^—1, the parameters a e have to be 
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large enough to ensure coerciveness of the bilinear form Bh on Xh (see, e.g., Lemma 2.1 of 
[22]). 

If u is a subset of fl, we denote by c a ^ (resp. C a ^) the minimal (resp. maximal) 
eigenvalue of the restriction of the matrix a to u. We further denote by cp^ u = inf xgti; (/i — 
^div/3). We recall that from our assumption, we have the implication 

c /3^,uj = =>■ div (3 = and /i = in u. 

For shortness we drop the index Q in these notations, namely we set c a = c Qi q, C a = C a ^ 
and cp^ = cp^ n . 

As our approximated scheme is a non conforming one (i.e. the solution does not belong 
to Hq(Q)), as usual we need to take into account the non conforming part of the error. 
The possibilities are twofold : either use an appropriate Helmholtz decomposition of the 
error (see Lemma 3.2 of [Hj, or Theorem 1 of [1] in 2D or Lemma 2.1 of [9]) or use the 
following comparison principle, proved in Theorem 6.3 of [T6j . 



Lemma 2.1 (Abstract upper error bound) Let u E Hq(VI) be a solution of (TJ)) and 
Uh G Xh be the solution of (TJ]) ; then 

\\W — Uh\\\< inf illlwft — sill (6) 

+ inf sup ( (/ — div t — j3 ■ Vs — /is, ip) 

te-ff(div,n) ^effi^ni^ii^i V 

- {aV h u h + t, Vcp) + ((/i - ^div/3)(s - u h ), cpj^ J. 

3 The a posteriori error analysis based on gradient 
recovery by averaging 

Error estimators can be constructed in many different ways as, for example, using residual 
type error estimators which measure locally the jump of the discrete flux |22j . A different 
method, based on equilibrated fluxes, consists in solving local Neumann boundary value 
problems [2] or in using Raviart-Thomas interpolant [TJ EJ HHJ [16]. Here, as an alternative 
we introduce a gradient recovery by averaging and define an error estimator based on a 
if (div )-conforming approximation of the broken gradient aVhUh- In comparison with |18j . 
we admit discontinuous diffusion coefficients and use a discontinuous Galerkin method. 

Inspired from [18] the conforming part of the estimator t]cf involves the difference 
between the broken gradient aWhUh and its smoothed version Guh, where Guh is for the 
moment any element in X\ 1 satisfying 

Gu h e H(div,Q) = {v G L 2 {nf : div v G L 2 (fi)}, (7) 
{Gu h )\\n 3 GH 1 (n j ),Vj = J. (8) 
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Hence conforming part of the estimator tjcf is defined by 



V 2 cf 



Yl ^cf,t, (9) 



TeT 

where the indicator i]cf,t is defined by 

Vcf,t = \\a~ 1/2 (aVwh - Gu h ) \\ T . 

For the nonconforming part of the error, we associate with Uh, its Oswald interpolation 
operator, namely the unique element Io s Uh G Sh defined in the following natural way (see 
Theorem 2.2 of [22]): to each node n of the mesh corresponding to the Lagrangian-type 
degrees of freedom of Sh, the value of Io s u h is the average of the values of Uh at this node 
n if it belongs to f2 (i.e., Io s Uh{n) = ^"'^ n ' > ) and is zero at this node if it belongs to 
T. Then the non conforming indicator t]nc,t is simply 

Vnc,t = \\\Io a u h - u h \\\ T = (||a 1/2 V(/ s^ - UhjWr + IK/" ~ ^div (3) 1/2 (I Qs u h - u h )\\^ . 
The non conforming part of the estimator is then 

V 2 NC=J2V 2 NC,T- ( 10 ) 
TeT 

Similarly by keeping in t)nc only the zero order term, we get a second non conforming 
indicator t] NC 2,t, which is simply 



Vnc2,t = 1 1 (// - T^div (3) 1/2 {I 0s Uh ~ u h ) 



with a global contribution 



VNC2 



^2v 2 nc2,t- (n) 



TGT 

We obviously notice that 

VNC2,T < Vnc,t- 

Similarly we introduce the estimator corresponding to jumps of Uh- 



with 



Vie 



h^lluhWl if e G £ D . 



As in [TS], we introduce some additional superconvergent security parts. In order to 
define them properly we recall that for a node x G J\f, we denote by A x the standard hat 
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function (defined as the unique element in Sh,i such that X x (y) = S x , y for all y G M), let 
u x be the patch associated with x, which is simply the support of X x and let h x be the 
diameter of u x (which is equivalent to the diameter Hk of any triangle K included into 
u x ). We now denote by r the element residual 

r = / + div (Gu h ) - (3 ■ VI 0s u h - fiI 0s Uh 

and for all x E M, we set 

r x = {L^xT 1 J w jK if x G Mint = M\Md, 
f x = o if x e M D - 



We further use a multilevel decomposition of S^i, namely we suppose that we start from 
a coarse grid 7o and that the successive triangulations are obtained by using the bisection 
method, see [HI [32]. This means that we obtain a finite sequence of nested triangulations 
%, £ = 0, • • • , L such that 71 = T ■ Denoting by Se the space 

S t = {veC(n)\v [ TeWi(T),TeTt}, 

then we have 

Se C Se + i and Sh,i = uf =0 S^ = Sl- 
Furthermore if we denote by Mi the nodes of the triangulation %, we have 

Af £ c W m . 

As usual for all z G Mt we denote by \t z the hat function associated with z, namely the 
unique element in Se such that 

Afe(z') = b zz Nz' EM t . 

For all i > 1 we finally set 

Me = (Me \ Me-i) U {z e Me-i ■ A fe ^ X e ^ lz }, 

and Mo = Mo- It should be noticed (see for instance [18]) that to each z e Me, the 
corresponding hat function Xe z does not belong to Si-i. 
Now we define p and 7 by 
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f = E E ^ 2 =. 
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where 



Fi "'a; / I' 'a:| ^xi 



llz = \(R,^iz)\, 

R being the residual defined by 

(R, (p) = (Gu h -Vcp + ■ VI 0s u h ip + pl 0s u h ip - f(p) ,V(p G H 1 (Q) . 
Jn 

3.1 Upper bound 

Theorem 3.1 Let u G Hq(Q) be a solution of problem and let be its discontinuous 
Galerkin approximation, i.e. Uh G Xh solution of O. T/ien i/iere exists C > swc/i i/iai 

w/illl < ?7cf + Wc + Vnc-2 + C(p + j), (12) 

and consequently 

\\u - u h \\ DGjh < ricF + Vnc + Vnc2 + VJ + G{p + j). (13) 
Proof: We first use the estimate (jSJ) with s = Io s u h and t = —Gu h to get 

WW - U h \\\ < {\\\u h - i 0a u h \\\ 

+ sup ( (/ + div Gu h - (3 ■ VI Qs u h - pl 0s u h , tp) 



<peHi(Cl):\\\<p 

aV h u h - Gu h ,Vip) + ((// - ^div 0)(I Os u h - u h ),cp) 



2 

By Cauchy-Schwarz's inequality we directly obtain 



\u - u h \\\ < tjcf + Vnc + VNC2 + sup \{R,f}\- (14) 

v eHi(ay.\\\<p\\\=i 



Using the arguments from Theorem 4.1 of [18J, we have 

\(R,ip)\<C{p + i)\W 1/2 Vip\\. (15) 
These two estimates lead to the conclusion. ■ 

Remark 3.2 Let us notice that under a superconvergence property of | \a~ 1 ^ 2 (Guf l — 
dVu)\\, p and 7 will be proved to be negligible quantities (see Theorem 13.81 below), so 
that the error is, in this case, asymptotically bounded above by the estimator without any 
multiplicative constant. This superconvergence property is observed in most of practical 
cases, as for example in our numerical tests (see section [5]). Moreover, theoretical results 
for different continuous finite element methods on structured and unstructured meshes 
have been established (see for example [T9l EU [39]), but, to our knowledge, not yet for 
discontinuous methods for reaction-diffusion-convection problems on unstructured multi- 
dimensional meshes. ■ 
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3.2 Asymptotic nondeterioration of the smoothed gradient 

This subsection is mainly devoted to the estimate of the error between the smoothed 
gradient and the exact solution, where we show that it is bounded by the local error in the 
-DG-norm up to an oscillating term. 

We first start with the estimate of some element residuals. On each element T, let us 
introduce 

R T = f + div (a T Vu h ) - (3 ■ Vu h - fiu h . 
Lemma 3.3 For all T E T , one has 

hrWRrh < (^VM^I^ + ^^^llk-^lllT + ^IIZ-nr/llr, (16) 

C a,T C /3,^,T 

h T \\r\\ T < Cl^a^l^aVu-Gu^T + hr^^ 

C a,T C P,fJ,,T 

+ M/-nT/lh 

where Hrf is the L 2 (T)- orthogonal projection of f onto F (T). 

Proof: We start by the first estimate. Its proof is quite standard, we give it for the sake 
of completeness. Denoting by 

R' T = Ii T f + div (a T Vu h ) - (3 ■ Vu h - /iu h , 

we trivially have 

R T = f - U T f + R' T . 

Hence it remains to estimate /irll-Rrllr- F° r that purpose, we introduce the standard 
bubble function (see [35J) and set Wj> = b^R! T . Then by a standard inverse inequality, 
we have 

ll-^rllr ~ / Rt^t 

JT 

< J (U T f + div (a T Vu h ) - /3 ■ Vu h - fiu h )w T 

< / (n T J - f)w T + (f + div (a T Vu h ) - f3 ■ Vu h - /2u h )w T . 

JT JT 

Now using dl]), we may write 

\\ r 't\\t ~ J ( n T/ - f)wr + J (div (a T V (u h - u)) -(3- V{u h -u) - fj,(u h - uj^w T 
and by Green's formula, we obtain 
\\R't\\t< J \n T f-f)w T - J (a T V(u h -u)-VwT + (/3-V{u h -u)+v{u h -u))w T y (18) 
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By Cauchy-Schwarz's inequality we obtain 

WKWt < l|n T /-/||TlkT||T + ||a^ /2 vK-«)|| T c'^ / T 2 ||Vw T || T 

+ ( ||/3||oo,TC^ /2 ||ar 2v ( M ^ ~ u )\\t + ll^l|oo,TC^||(/i- \d\v f3) 1/2 (u h - u)\\ T J ||wr||r. 



The inverse inequalities [35] 

\\w t \\t <\\R' t \\t, (19) 
\\Vw T \\ T <h^\\R' T \\ T , (20) 

lead to IjlSj) . 

We proceed similarly for the estimate (1TT1) . namely for any T G T we set 

= n T / + div - /3 ■ VI 0s u h - ju/os^/i, 

and therefore 

r = f — HtJ + on T. 

Hence it remains to estimate ^tII^tIIt- As before we set wt = bj^vip. Then by a standard 
inverse inequality, we have 



-LWl < 



\\ r T\\T ,b / r T w T 

Jt 

< / (II r / - f)w T + / (/ + div (Citft) - f3 ■ VI 0s u h - fiI 0s u h )w T . 
Jt Jt 

Using (JT]) we obtain 

H r TllT< J ( n r/ - f)w T + J (div (Gu h - a T Vu) - (3 ■ V(I 0s u h - u) - ^{I 0s u h - u]^w T 
and Green's formula yields 

ll r Tllr~ J ( n r/ - f)wr- J ({Gu h -a T Vu)-Vw T + (/3-V(Io s u h -u)+ n{I 0s u h -u))w T 



(21) 

Cauchy-Schwarz's inequality and the inverse inequalities ( TT9|) -( l20|) yield (jT7|) . ■ 

We go on with the edge residual. 
Lemma 3.4 For all e G Si nt , we set 

[aV h Uh -nj e = a T +Wu h \ T + ■ n T + + a T -Vu h \ T - ■ n T -, 
when uj e = T + U T~ . Then one has 

\\/3\\oo,T . HHIoo/Tn 



/#*||[aV h u fc .n]J e <^ (22) 
where osc(/» 2 = ^ Tcw h 2 T \\f - Ii T f\\ 2 T . 



C a,T C 0,n,T 
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Proof: Again the proof is quite standard, we introduce the edge bubble function b e and 
set w e = E(r e )b e , where E{r e ) is an extension of r e = [aVhUh • n] e to u e (see [35]). By a 
standard inverse inequality, we have 




By Green's formula, we obtain 

IMIe ~ / ( a T^(uh - u) ■ Vw e + div (a T V(u h - u))w e ) . 

Taking into account we get 

(a-rVtw/i — u) • Vw e + div (axV«/ l )io e + (f — (3 ■ Vu — jiu)w e ) 

(a T V(u h - u) ■ Vw e + (/ + div (a T Vu h ) - ■ Vu h - fiu h )w e 
{u h -u) + fi(u h - u))w e J 

(a T V(u h - u) ■ Vw e + R T w e + (/3 ■ V(u h - u) + /i(u h - u))w e ) . 

The conclusion follows from Cauchy-Schwarz's inequality, ( IT6l) and the next inverse in- 
equalities 

lke||T<^ /2 ||r e || e ,VTCu; e , 
l|Vw; e ||T<^ 1/2 ||r e || e ,VTca; e . 

■ 

To prove our nondeterioration result we also need the next estimate of the norm of 
the difference of a vector field v with an appropriate projection g(v) with the norm of its 
jumps in the interfaces. First of all for any vertex x of one ttj and belonging to more than 
one sub-domain, we introduce the following local notation: let fij, % = 1, • • • , n, n > 2, the 
sub-domains that have x as vertex. We further denote by the unit normal vector along 
the interface Jj between fij and fij+i (modulo n if x is inside the domain Q) and oriented 
from fli and Qi+i. Now we are able to state the following lemma, for a proof see Lemma 
3.3 of [T3]. 

Lemma 3.5 Assume that x is a vertex of one Qj and belonging to more than one sub- 
domain, and use the notations introduced above. Then there exists a positive constant C 
that depends only on the geometrical situation of the Vti 's near x such that for all £ IR 2 , 
i — 1, • • ■ , n, there exist vectors g{v)^ 1 ' £ R 2 , i — 1, • • ■ , n satisfying 

(g(v)^ - g ( v )M) . m = 0, Vi = 1, • • • ,n, (23) 



e ~ ^2 I 

~ E L 
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and such that the following estimate holds 



5> (0 -iK«) W l<C£l[ v '"I*!' (24) 

i=l i=l 

where here \ ■ \ means the Euclidean norm and \v ■ n\ { means the jump of the normal 
derivative of v along the interface lii 

\vn\ = -v®) -m^i = 1, 

Using the above lemma, we are now able to prove the asymptotic nondeterioration of 
the smoothed gradient if the following choice for Guh is made (we refer to p2] for a similar 
construction): We distinguish the following different possibilities for x E M . 

1) First for all vertex x of the mesh (i.e. vertex of at least one triangle) such that x is 
inside one Qj, we set 

(Gu h ) lQj (x) = y2\T\a T Vu hlT (x). (25) 
\uj T \ * — ' 

1 x ' x€T 

2) Second if x belongs to the boundary of Q and to the boundary of only one Qj (hence it 
does not belong to the boundary of another Q k ), we define (Guh) \nj(x) as before. 

3) If x belongs to an interface between two different sub-domain Qj and Q k but is not a 
vertex of these sub-domains, then we denote by n^ k the unit normal vector pointing from 
Qj to Qk and set tj >k the unit orthogonal vector of rij >k so that (nj tk ,tj jk ) is a direct basis 
of M 2 ; in that case we set 

{Gu h )\ aj (x) ■ n j>k 
(Gu h ) ln .(x) ■ t jtk 
(Gu h )\ Uk (x) ■ t j)k 

4) Finally if x is a vertex of at least two sub-domains Qj, for the sake of simplicity we 
suppose that each triangle T having x as vertex is included into one Qj, and we take 

{Gu h ) ]Qj {x) = g{v)®Vj e J x , (29) 

where J x = {j E {1, • • • , J} : x E f2j}, g(v)^' were defined in the previous Lemma I3T51 
with here v given by v = {ajVuh\T{x))j & j x . 
With these choices, we take 

(Guh)^ = ^ ( Gu h)\n 3 (x)X x ,yj = 1, • • ■ , J, (30) 



(Gu h )\Q k (x) ■ n jtk = - — r V" \T\a T Vu h \ T {x) ■ n jtk , (26) 

1 1 x&T 

I * | J2 \T\a T Vu hlT (x)-t^ (27) 



rcn,:xGT 



\UJ X fl Q k 



^— Yl \T\a T Vu hlT (x)-t j>k . (28) 



Tcn k :x€T 
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where (Guh)\n.(x) was defined before. 

The main point is that by construction Gu^ satisfies the requirements ([7]) and (jHJ) but 
moreover the next asymptotic nondeterioration result holds. In order to track the constants 
with respect to the diffusion coefficient we first prove the following technical lemma. 

Lemma 3.6 If{n,t} is an arbitrary orthonormal basis ofM 2 , then for allT E T, one has 

\a T v\ < K(a T )(\(a T v) ■ n\ + 2C a , T \v ■ t\), e M 2 , (31) 

where we recall that k(ot) is the condition number of the matrix ay defined by 

k(ot) = C a>T c^ T . 

Proof: For shortness we drop the index T. 

In a first step we take !»6l 2 such that w ■ t = 0. Then in that case we can write 

w = {w ■ n)n, 

and therefore 

aw = (w ■ n)an, aw • n — (w • n){an ■ n). 
By direct calculations, we obtain 

\aw\ < C a \w ■ n\, \w ■ n\ < c~ l \aw ■ n\, 

that leads to 

\aw\ < K(a)\aw ■ n\. (32) 

For an arbitrary v E M 2 we set w = v — (v ■ t)t that, by construction, satisfies w -t = 0. 
Hence by (|3"2~|) we get 

\av\ < \aw\ + \v ■ t\\at\ 

< K,(a)\aw • n\ + C a \v • t\ 

But by the definition of w, one has 

aw ■ n — av ■ n— (v ■ t)at ■ n, 

and consequently 

\aw ■ n\ < \av ■ n\ + C a \v ■ t\. 

The last inequalities leads to f )3T]) reminding that n(a) > 1. ■ 

Thanks to this lemma we can prove an asymptotic nondeterioration result that is similar 
to Theorem 3.4 of fT3], where we here give the explicit dependence of the constants with 
respect to the coefficients a, (3 and \x. 
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Theorem 3.7 If £ < 2, then for each element T G T the following estimate holds 

\\a- 1/2 (Gu h - a T Vu)\\ T < c^ 2 C a , T K(a T ) ^ C /2 |IMI| e + (33) 

e€£i„t:eCtUT 

1 i -1/2 / \ ( rH" 2 - , l /ll^lloo.T' . UmIIoo.T' x\ \ I,, mi 

1 + c a ^ K{a T ) max [CJ T , + h T { — + — ^ ) J I - w ft ||| WT 

— 1/2 

+c a T 7 «(ar) osc(/, uj t ), 

where ut denotes the patch consisting of all the triangles of T having a nonempty inter- 
section with T. 

Proof: By the triangle inequality we may write 

/2 — -1/2 1/2 

\\a T (Gu h - a T Vu)\\ T < \\a T {Gu h - a T Vu h )\\ T + \\a T (a T Vu h - a T Vu)\\ T 

— 1/2 

< \\a T (Gu h - a T Vu h )\\ T + \\\u- u h \\\ T . 

Therefore it remains to estimate the first term of this right-hand side. For that purpose, 
since T C Clj for a unique j G {1, • • • , J}, we may write having in mind the assumption 
I < 2 : 

(Gu h - a T Vu h )\ T = ^{(G« ft )^.(a;) - a j Vu h \ T (x)}\ x . 

As < X x < 1, and since the triangulation is regular, we get 

1/2 (Gu h - a T Vu h )\\ T < c^ 2 \(Gu h )\ aj {x) - aj Vu h \ T {x)\h T . (34) 



\ a T 



We are then reduced to estimate the factor \{Guh)\o,j{.x) — ajVuh\T(x)\ for all nodes x of 
T. For that purpose, we distinguish four different cases: 

1) If x G Qj, then we use an argument similar to the one from Proposition 4.2 of [18] 
adapted to the DG situation. By the definition of Guh, we have 



(Gu h )\ a (x) = —- Y] \T'\ajVu h \ T/ \ 
\Ux\ _r — 



x) 



because in this case all T' C u x are included into Qj. As a consequence, we obtain 
{Guhjfoix) - ajVu h \ T (x) = | — r \T'\aj(Vu h \T'{x) - Vu h \ T (x)), 

and therefore 

\(Gu h )\n 3 (x) - ajVu h \ T (x)\ < ^ ^(Vi^T'Or) - Vu h \ T (x))\. 

T'Cuix 
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For each T' C u x , there exists a path of triangles of cu x , written Ti, % — 0, • • • , n such that 

T = T,T n = T',T t ^T J ,\/t^j, 
Ti fl T i+ i is an common edge Vi = 1, • • • , n — 1. 

Hence by the triangle inequality we can estimate 

n-l 

\aj(Vu h \ T '{x) - Vu h \ T (x))\ < ^2 \ a j(^ u h\T l+1 (x) - Vu h \ Ti {x))\. 

Now for each term, since cij is symmetric and positive definite, using Lemma 13.61 above we 
have 

|aj(Vu h |T (+1 (a:) - Vu h \ Ti (x))\ < <ajOI{ a i(Vufc|T i+ i(aO - Vu^a))} • rn\ 

+2K(a j )C aj \(Vu h \T t+1 {x) - Vu h \ Tt (x)) ■ Ul 

where rii is one fixed unit normal vector along the edge Tj D Tj + i and ti is one fixed unit 
tangent vector along this edge. All together we have shown that 

\(Gu h )\ nj (x) - ajVu h \ T (x)\h T < h T K(a T ) ^ {\\aVu h (x) ■ n\ e \ + C a ,r\ IV^(x) • tj e \}. 

e&£i nt :x£e 

Using a norm equivalence and an inverse inequality we obtain 
\(Gu h ) ]Qj (x) - aj Vu hlT (x)\hr < k(ot) ^ {/# 2 ||[aV fc u A ■ n]|J c + C^/iJ^H 

(35) 

2) If the node x belongs to the boundary of Q and to the boundary of a unique ttj, since 
(Guhj^jix) is defined as in the first case, the above arguments lead to (1351) . 

3) If x belongs to an interface between two subdomains and is not a vertex of them, then 
it is not difficult to show that 

\{Gu h )\ aj (x) - ajVu h \T(x)\hT < h T ^ |[aV h u fc (x) • n]J (36) 

holds (due to the regularity of the mesh), and consequently (1331) is still valid. 

4) Finally if x is a vertex of different sub-domains Qj, then Lemma 13.51 yields ( 13T)|) and 
therefore as before we conclude that (13"5]) holds. 

Summarizing the different cases, by (134|) and (1331) . we have 

||a- 1/2 (G Mh -a r V^)|| r <c^ / \(a T )^ £ {^ /2 || [aV h u h • n] Je+Oa^llKllle}- 

xGT eG£i 

(37) 

The first term of this right hand side is the standard edge residuals that were estimated 
in Lemma 13.41 while the second term in ( 1371) is part of the DG-norm and is then kept. 
Therefore the estimate (p2l) in (1371) leads to ([33]). ■ 
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3.3 Lower bound 

In the spirit of subsection 13.11 (see Remark 13. 2p . we first provide lower bounds where the 
error between the gradient of the exact solution and its smoothed gradient is involved in 
the right-hand side. In a second step using the results of the previous section, we give 
lower bounds with only the DG-norm of the error. 

First using the same arguments than in Proposition 4.1 of [18], we have 

Theorem 3.8 For all T G T, x G M and t > 0, z G Mi, we have 

Vcf,t < \\a](\(u h - u)\\ T + \\a~ l/2 (Gu h - aVu)\\ T , (38) 

c-X 2 p x < ^\\a-^(aVu-Gu h )\l x (39) 

C a ,oj x 



+ h x I + 1/2 ye )\W U ~ IosUh\\ k + c a J x osc(/, U x ), 

v C a ,ui x Ca,u~C- ' 



It. < ^IM- 172 ^^ - + ^(^^- + 111^ - ^os^lll^- (40) 

Ca > z C /3,li,z 

Proof: The first estimate is a simple consequence of the triangle inequality. For the second 
one, we notice that 

Px — h x 1 1 v 1 1 UJx . 

Hence the estimate (l39l follows from (fl7]l . 

For the third estimate by the definition of ji z , ([T]) and Green's formula, we have 



llz 



I {{Gu h - aVu) ■ VA fe + (J} ■ V{I 0s Uh -u) + n{I 0s Uh ~ u))X iz ) . 
Jn 

Hence the conclusion follows from Cauchy-Schwarz's inequality and using the estimates 

< fits, llVAfell < 1. 



For the non conforming part of the estimator, we make use of Lemma 3.5 of [16] (see 
also Theorem 2.2 of [22]) to directly obtain the 

Theorem 3.9 Let the assumptions of Theorem \3.1\ be satisfied. For each element T G T 
the following estimate holds 

Wa > T<« / T 2 + /i T |K/.-^div/3)^|| 00)T ) K 1/2 \\bhj\\ e . (41) 
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Corollary 3.10 Let the assumptions of Theorem \3.1\ be satisfied. For each element T G T 
the following estimate holds 

Ilk - IosUhWW < (1 + C l J} + MIG" " 2 div / 3 ) 1/2 H-^)ll M - M\dgw (42) 

where 




Remark 3.11 From the estimates (140)) and (]42p . we can say that the quantities 

CaJJxpx, P and are superconvergent if ||a _1 / 2 (aVw — Crtt^) || and osc(/, u^) are. Obviously 
this superconvergence properties are lost when the diffusion matrix a tends to 0, i.e., it is 
not robust with respect to convection dominance. A similar phenomenon also holds for 7, 
see Remark [3. 151 I 



A direct consequence of these three Theorems is the next local lower bound: 

Theorem 3.12 Let the assumptions of Theorems \3.1\ and 3/7_ be satisfied. For each ele- 
ment T e T the following estimate holds 



VCF,T + VNC;T + Vj,t + ^{c a ^£ Px + 7x0 < «i,t||« - u h \\ DGi z T 

x£T 

+c«i 2 T ( 1 + C?) IK 1/2 («Vu - Gu h ) |U + E osc ^ ^ 

where 7^ = ^Lx recalling that L is such that Ml = N, Cot is the larger patch defined by 
Cj t = Ut'cui t ^t' and 

= (1 + CaS T + M<A* - \div^\\^ T ){l + h T (l + C-J/ 2 ) (Ife + )) }. 

Using the nondeterioration result from the previous subsection, we get a lower bound 
with the .DG-norm of the error at any event. 

Corollary 3.13 Let the assumptions of Theorems \3. 1\ and 3. 7 be satisfied. For each ele- 
ment T G T the following estimate holds 

Vcf,t + Vnc,t + Vj,t + ^2(c~^p x + 7,.) < - u h \\ DG ^ T + c~^k(o Wt ) 2 osc(/, wy), 

where 



*1,T = «!, T + C^l + C-V 2 ) (l + 



2 /• \2 



+ C a,u T K \ a uT)\ Lj a,u T + M Y/2 1 V2 >> I' 

c c 
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Combining the arguments of Proposition 4.3 of [18J and the above ones, one can obtain 
a global lower bound: 

Theorem 3.14 Let the assumptions of Theorems \3.1\ and be satisfied. Then the fol- 
lowing global lower bound holds 

Vcf + Vnc + Vj + P + l ^ ^2 1|^ - u h \\ D G,h + osci(/,fi), 

where 

«2 = «3 

"V2^ _^„_^ , „-l/2^r„_> max ^V2 , ^_JI/3||oo,T' , ||/X|| 

T'Cw T 



f, , f -1/2^ / \ . -1/2 , \ ^1/2 , L / P oo,T' . AMoo,T'x\ \ 1 

+ + max cj, C a)T K{a T ) + c a + K{a T ) max [C^ T , + /i T ( — ^— + — ) J J, 

V " T C a,T' C AjU,T' / 

«3 = (1 + C Q 1/2 ) (l + max/^(Ml^ + Nkg, ) + c-fdl/^IU + Hdiv^lU + ||^|U)) , 



I v 1/2 / x 1/2 

OSC 1 (/, Q) = I ^ C"^ 0SC(/, W z ) 2 I + K l I XI C a,T K ( fl T) 2 OSc(/, W T ) 

VzeAf / \TeT 



Proof: As in Proposition 4.3 of [18], we have 

with x £ -f^o(^) defined in the proof of Proposition 4.3 of [18]. As before using (j2J), we 
have 

x) = / ((Gufc - aVu) ■ V X + (P ■ V(I 0s u h - u) + n{I 0a u h - u))x) (43) 
Jn 

and by Green's formula, we get 

(R,X)= ((Gu h -aVu)-Vx+(Io 8 Uh-u)div(f3x) + Klo B Uh-u)x)- (44) 
Jn 

By Cauchy-Schwarz's inequality we get 

(R,X) <C a 1 / 2 ||a- 1 / 2 (Gn,-aVn)||||Vx|| + (||/3|U + Hdiv/^IU + H/^IU) ||« - / 0s ^||||x||i,n. 
Hence by Poincare's inequality, we deduce that 

f < {c l J 2 \\a~ l/2 (Gu h - aVu)\\ + (j^lU + \\diy ^ + |MU)||« - ios«h||) ||V X ||. 
Since Proposition 4.3 of [18] shows that ||Vx|| ^ 7; we conclude that 

7 < C a 1/2 ||a^ 1/2 (G^ - aVu)\\ + (ll/^ + Hdiv^U + |MU)||« - W^ll- (45) 

This estimate and Theorems 13.71 13.81 and 13.91 lead to the conclusion. ■ 
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Remark 3.15 From the estimate (1451) we can say that the quantity 7 is superconvergent if 
||a _1 / 2 (aVM— Guh)\\ and Io a Uh\\ are superconvergent. For the second term \\u— Io s Uh\\, 
using a triangle inequality, we have 

\\U - JosUftH < \\U - U h \\ + \\U h - IosUh\\, 

and Lemma 3.5 of [16] yields (see Theorem 13.91) 

IK - Io s Uh\\ < h\\(ii - ^div f3) 1/2 \\ 00 \\u - u h \\ DG , 

Hence \\uh — Ios u h\\ is a superconvergent quantity. On the other hand using an Aubin- 
Nitsche trick (see for instance [3]), if 9 = 1, then we have 

\\u - u h \\ < Ch a \\u - u h \\ DG , 

for some o > and C depends on the matrix a, f3 and fi (and may blow up as a — > 0). 
Hence, if 9 = 1, \\u — Uh\\ is is a superconvergent quantity and therefore \\u — Io s Uh\\ as 
well. As in Remark 13. 11[ this property is not valid in the convective dominant case. ■ 



According to Theorem 13.141 we see that our lower bound is not robust with respect to 
convection/reaction dominance. Indeed two factors (namely the factors c^JI/^Hoo^ and 

c a ,i{ 2 1 1 /x 1 1 oo,lu x in the definition of k 3 ) blow up as a goes to zero. The two factors come from 
two different terms in the proof of the lower bound. Indeed the first one appears when 
one wants to estimate a convective derivative (see Lemma 13.31 for instance) and is usually 
eliminated by adding to the norm of the error an additional dual norm (see [531 ES]), this 
technique will be adapted to our problem in the next section. We further see that the same 
factors appear in the estimate of the error between Gun and aVw, hence a simple idea to 
eliminate these factors is to add the term | |a~ 1//2 (G , M/ i — dVu) \\ to the error (even if in 
practice this error is quite often superconvergent). These two ideas are developed below. 

Note also that in Theorem 13.141 the lower bound is not robust with respect to the 
local anisotropy of the diffusion matrix (via the constants n\ and K2 that blow up if the 
condition number n(cij) blows up for at least one j). But this is a normal phenomenon 
that also appears in former works, see for instance Theorem 3.2 of [16], Theorem 4.4 of 
[37] or Theorem 4.2 of [38]. 

4 A robust a posteriori estimate 

According to the previous papers [331 EH] we need to use an appropriate dual norm in 
order to estimate convective derivatives. The price to pay is that the constant in the upper 
bound will be no more 1, but remains nevertheless explicit. Note further that theoretically 
we lose the robust superconvergence property of p (or a modification of it, see below) and 
of 7. 
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We start with the following definition: for h £ L 2 (Q), let us denote by its norm 

as an element of the dual of Hq(Q) (equipped with the norm ||| • |||), namely 

= sup ffr-TTT- 

veHi(n)\{o} IIMII 

Then the proof of Lemma 3.1 of [36] yields the following result. 
Lemma 4.1 If 



M = max{l, sup f }, 

n (/i-sdiv/?) 



i/ien we have 



(aVv ■ Vw + fxvw) 



< M\\\v\\ \ \\\w\\\,Vv,we H%(n), 



and 

|||/3-V™|||*<(2 + M) sup ^Vwe^n). 

«GHi(n)\{0} 1 1 Fill 



\w\ 



For robustness reasons, we further introduce 



a x = mm{c^ 2 ,c a ^h x },\/x e AT, 
a T = xmn{c^ 2 ,c-^ 2 h T },WeT, 



and 



J2 a * \ r ~ r x\ 2 K- 



Theorem 4.2 Let u £ i?o(^) ^ e a solution of problem (TJ|) and Zei «/i &e its discontinuous 
Galerkin approximation, i.e. Uh £ Xh solution of (TJ]). Then there exists C > sttc/i t/iat 

|||M-u h ||| + HI/9- V(u-IosOI||* < + (3 + M) (r/ C F + Wc + C(p + 7)), (46) 
and consequently 

\\u- u h \\ DG>h + 1 1 1/3 • V(« - /os^ft) HI* < + ^2 + (3 + M) (r/ CF + + C(p + 7) ) . (47) 
Proof: Applying Lemma [4.11 to w = u — Io s u h: we see that 

\\\P-V(u-I 0s u h )\\U<(2 + M) sup i?C "~, / ^^ ) . 

"eff„ 1 (n)\{0} 1 1 Fill 

But according to (j2j) and the definition of i? we have 

B(u-I 0s u h ,v) = {a- 1/2 (Gu h -aV h u h ),a- 1/2 Vv) + {a 1/2 {V h u h -VI 0s u h ),a- 1/2 Vv)-(R,v). 
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Hence by Cauchy-Schwarz's inequality we deduce that 

B(u - I 0s u h ,v) < {ricF + TMrcOIIMH + \{R,v)\. 
Similarly to Theorem 4.1 of |18j . we can show that 

mv)\<c{p+i)\\\<p\\\, (48) 

due to the easily checked estimate 

\\<p-(f K) I <pK\\ < a^HMIL- 



The estimate ( 148]) then yields 

B(u - I 0s u h , v) < (rjcF + Vnc + C(p + 
and therefore 

\\\P-V(u- I 0s u h )\\U < (2 + M)( VC F + VNC + C(p + *f)). (49) 
In the same manner, the estimate (1481) allows to show that 

\\\u - u h \\\ < r] C F + Vnc + VNC2 + C(p + 7). 
This estimate and fl49l) lead to fl46l). ■ 



As suggested before we now add the error between the gradient of u and its smoothed 
gradient in order to have a full robust lower bound. 

Corollary 4.3 Under the assumptions of Theorem \4-2\ there exists C > such that 

|||u-Ufc||| + \\a~ 1/2 (Gu h -aVu)|| + ■ V(u - Io s u h )\\\* < VcF + 2r] NC2 

+2(3 + M) ( Vc f + Vnc + C(p + 7)) , 

and consequently 

\\u - u h \\ DG ,h + \\a~ 1/2 (Gu h - aVu)\\ + \\\/3 ■ V(n - J s^)|||* < Vj + Vcf + %Vnc2 

+2(3 + M) {vcf + Vnc + C(p + 7)) . 

The lower bound requires to revisit some results from subsection 13.31 For that purpose 
and only for the sake of simplicity we require that the diffusion matrix a is constant on 
the whole domain. Let us now revisit Lemma 13.31 

Lemma 4.4 Assume that the diffusion matrix a is constant. Then the next estimate holds 



p < ||a- 1/2 (G^ - aVn)|| + \\\/3 ■ V(« - /osOIH* + ^7? \\u - I 0s u h \\ + £, (50) 
where we have set 



e = Y, a T\\f-^f\\ T - 



TeT 
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Proof: First we notice that 



xeM TeT 



Hence it remains to estimate 



. 2 1 1 2 

C 

TeT 



Oirp 1 1 1 1 T" 1 • 



For that purpose, we start as in the proof of Lemma 13.31 Namely for all elements T, the 
estimate (l2Tj) holds. Hence multiplying this estimate by a^, summing on T and setting 
w = J^TeT a T w T (with wt defined in the proof of Lemma I3~3l) . we have 

E^KIIt < E«T / (^Tf-f)w T 
TeT TeT ^ T 

- J (Gu h - aVu) ■ Vw + ((3 ■ V(I 0s u h - u)w + f^{Io s Uh - u)w^j . 
By Cauchy-Schwarz's inequality and the definition of the norm ||| ■ HI*, we get 
E^IKIIt < e(J2 a T\M\T) + (h~ 1/2 (Gu h -aVu)\\ 



TeT TeT 



j —Til' TllT' 
TeT 



+ HI/3 • V(lt - ios«/i)|||* + ||tt - ios«h||)IIIH 

By the estimate (4.15) of [36], we see that 

iNir<xXn^^ 

while we directly check that 

TeT TeT 

These three estimates yield 



2 ||2 



(5> t |KIIt) £ l|a- 1 ^fc-aVu)|| + ||| j 9-V(«-J 0i « fc )|||. + ^||u-/o.«fc||+e, 

TeT c P,ij, 

and we obtain (1501) by the definition of r^. ■ 

As in Theorem 13.144 the above results allow to obtain the following robust global lower 
bound. 
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Theorem 4.5 Let the assumptions of Theorems\3 . 1\ and\3.7\ be satisfied. Assume that the 



diffusion matrix a is constant. Then the following global lower bound holds 

VCF + VNC + VJ + P + 7 ^ K A \\\u - U h \\ \ + K 5 T]j 

+«6|||/3 • V(u - / 0s ^)lll* + h~ 1/2 (Gu h - aVu)\\ + «(o)(l + C l J 2 )£, 

where 

K A = M(l + c^), 



l + c^ + iK/i-^div/?) 1 / 2 ! 



\x\ 



Proof: By using the identity P3l) of the proof of Theorem 13.141 we get 

f = (R,x) 

< C^Wa-^iGuu - aV«)||||Vx|| + IMUI« - IoM\\\x\\ + V(u- I s0lll* 
But Poincare's inequality yield 

iiixiir^^+ii/i-^div/jn^iivxii 2 . 

These two estimates and the fact ||Vx|| < 7 (see Proposition 4.3 of [18]) yield 
7 < Cy 2 ||a- 1/2 (G'« ft -aVw)|| + |^|| 00 ||it-/o s w/ l || 
+ (C a 1/2 + II (A* - ^iv /?) 1/2 |U) 1 1 1/3 • V(u - 7o.« fc ) I Ik 

This estimate, (1351) . (150]) and ( )4T|) lead to the conclusion. 
Remark 4.6 The factor M in the lower bound is also present in Theorem 4.1 of [36] . 



5 Numerical results 

Here we illustrate and validate our theoretical results by some computational experiments. 
5.1 The homogeneous case 

We consider the domain Q = {0 < x, y < 1}, the reaction coefficient \i — 1, the velocity 
field j3 = (1,0) T and the isotropic homogeneous diffusion tensor a = eX where X is the 
identity matrix. Here, as in [TB], we take e = IE — 02 and e = IE — 04. The source term 
/ is chosen accordingly so that u = \x(x — l)y(y — 1)(1 — tanh(10 — 20x)) is the exact 
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Figure 1: The exact solution u = —x{% — l)y(y — 1)(1 — tanh(10 — 20a;)). 



solution (see Figured]). 

Results are presented for uniformly refined meshes. In Tables Q] and El N stands for 
the number of mesh elements, rj = t\ C f + Vnc + VNC2 + VJj an d Eff = v/\\ u ~ u h\ \dg,h is 
the effect ivity index. CV error and CV recov are respectively the convergence order in \Jl/N 
of the error ||u — Uh\\DG,h and of || a -1 / 2 ((?«/, — aVu)||, from one line of the table to the 
following. First of all, it can be observed for e = IE — 02 that the error \\u — Uh\\DG,h 
converges towards zero at order one and that, in the same time, the super convergence 
property of | \a~ 1 ^ 2 (Guh — aVu)|| is observed. Moreover, as expected by Theorems 13.11 
and 13.141 when p and 7 are superconvergent terms, the proposed estimator is reliable and 
efficient since the effectivity index remains constant during the refinement process (around 
2.00). For e = IE — 04, the same kind of behaviour occurs. Let us nevertheless note that 
the convergence rate CV error is astonishment high. In fact, the contribution of the jump 

1 /2 

term (J2 e e£ ^e 1 II l u h} arising in \\u — u^Wdg^ is predominant and goes fast towards 
zero, while the contribution of — Uh\\\ is smaller but converges at order one once the 
mesh is fine enough. 

5.2 The singular case 

We consider here the domain Q = { — 1 < x, y < 1}, which is decomposed into 4 sub- 
domains fy, i = 1, 4, with fix = (0, 1) x (0, 1), Q 2 = (-1, 0) x (0, 1), fi 3 = (-1, 0) x (-1, 0) 
and f2 4 = (0, 1) x (—1, 0). Like in subsection 15. 11 the reaction coefficient is \i — 1, the veloc- 
ity field f3 = (1,0) T , but the isotropic diffusion tensor is this time no more homogeneous. 
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N 


\\ u — u h\\DG,h 


rv 

^ ' error 


r) 


Eff 


\\a-y 2 {Gu h -aVu)\\ 


rv 

^ ' recov 


512 


2.22E-02 




4.47E-02 


2.01 


4.29E-03 




2048 


1.13E-02 


0.98 


2.29E-02 


2.02 


1.61E-03 


1.41 


8192 


5.58E-03 


1.02 


1.13E-02 


2.02 


4.87E-04 


1.73 


32768 


2.77E-03 


1.04 


5.57E-03 


2.01 


1.37E-04 


1.83 


131072 


1.38E-03 


1.01 


2.77E-03 


2.01 


3.96E-05 


1.79 



Table 1: Homogeneous case, e = IE — 02, 7 = 250, 7 aj6 = e. 



N 


\\u - u h \\ DGjh 


rv 

^ v error 


V 


Eff 


\\a- l > 2 {Gu h -aVu)\\ 


rv 

^ v recov 


512 


2.60E+00 




2.72E+00 


1.05 


3.52E-02 




2048 


1.46E+00 


0.82 


8.28E-01 


1.06 


2.17E-02 


0.69 


8192 


5.38E-01 


1.44 


5.81E-01 


1.09 


7.97E-03 


1.44 


32768 


1.84E-01 


1.54 


2.07E-01 


1.12 


2.40E-03 


1.73 


131072 


6.83E-02 


1.43 


8.06E-02 


1.18 


7.14E-04 


1.75 



Table 2: Homogeneous case, e = IE — 04, 7 = 250, 7 aje = e. 



It is denned by a\Q. = £jX, with e 2 = £4 = I and E\ = e s = C > 1 to be specified. 

Using the usual polar coordinates (r, 9) centered at (0,0), the exact solution is chosen 
to be equal to the singular function u(x,y) = r](r) v(x,y) with v(x,y) = r a <f)(8), where 
a G (0, 1) and are chosen such that v is harmonic on each sub-domain % — 1, ..4, and 
satisfies the jump conditions : 

[v] = and [aVf.n] = 

on the interfaces. The function 77 is a C 1 [0, 1] truncation function used to ensure homoge- 
neous Dirichlet boudary condition on the boundary. Namely, we chose : 

( 1 for < r < 1/3, 

rj( r ) = < (r - 2/3) 2 (54r - 9) for 1/3 < r < 2/3, 
[0 for r > 2/3. 

It is easy to see (see for instance [11]) that a is the root of the transcendental equation 

a it 1 



tan 4 

and since a < 1, this solution has a singular behavior around the point (0, 0). Consequently, 
a local mesh-refinement strategy is used, based on the estimator t]t = tjcf,t + Vnc,t + Vj,t 
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and the marking process 

r] T > 0.75 max 777V, 

with a standard refinement procedure associated with a limitation on the minimal angle. 

Like in [T3], we choose first C = 5. Figure [2] shows some of the meshes obtained during 
the local refinement process. Moreover, Table [3] displays the corresponding quantitative 
results. It can be observed that the error goes towards zero as theoretically expected, and 
that the effectivity index always remains almost constant, which is quite satisfactory and 
comparable with results from [HI US] as well as those of the previous test in subsection 15.11 
The superconvergence property of | |a~ 1 / 2 (Gu/ l — aVu)|| is once again observed, and the 
mesh is automatically refined in the vicinity of the singularity as well as in the zone of the 
mesh where the gradient of 77 is the highest. 




Figure 2: Mesh levels 1, 3 and 5, singular solution for C — 5. 



N 


\\u - u h \\ DG;h 


rv 

^ v error 


77 


Eff 


\\a~^ 2 (Gu h - aVu)|| 


rv 

^ v recov 


128 


1.26E+00 




3.09E+00 


2.44 


1.14E+00 




468 


4.72E-01 


1.51 


1.25E+00 


2.65 


4.55E-01 


1.41 


2016 


2.14E-01 


1.08 


5.74E-01 


2.68 


1.36E-01 


1.65 


9068 


1.01E-01 


1.00 


2.81E-01 


2.76 


4.68E-02 


1.41 



Table 3: Homogeneous case, C = 5, 7 = 50, ^y^e = 1- 



Secondly, we consider a stronger singularity by choosing C = 100. Figure [3] shows 
some of the meshes obtained during the local refinement process, and Table H] displays the 
corresponding quantitative results. The results are very similar to the ones obtained for 
C = 5. Like in [13], the refinement process is faster around the interfaces (and the origin). 
The effectivity index slightly increases while remaining constant during all the refinement 
process, and the superconvergence property of | |a~ 1//2 (G ? M/ l — aVw)|| is observed. 
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Figure 3: Mesh levels 1, 3 and 15, singular solution for C = 100. 



N 


\\u - U h \\ DG ,h 


cv 

^ v error 


V 


Eff 


\\a'^ 2 {Gu h -aVu)\\ 


CV 

^ ' recov 


128 


1.56E+00 




6.34E+00 


4.05 


2.92E+00 




504 


7.25E-01 


1.12 


3.65E+00 


5.03 


7.99E-01 


1.89 


2208 


4.48E-01 


0.65 


1.95E+00 


4.35 


3.81E-01 


1.00 


7716 


2.44E-01 


0.97 


1.09E+00 


4.46 


1.60E-01 


1.38 



Table 4: Homogeneous case, C = 100, 7 = 500, 7 ai6 = 1. 



5.3 The boundary- layer case 

We now consider the domain Q = {0 < x, y < 1}, with the reaction coefficient /1 = 
and the velocity field f3 = (1,0) T . The homogeneous isotropic diffusion tensor is de- 
fined by a = el, with e = IE — 02. The source term / is chosen accordingly so that 
u = 10y(l — y)x(e~ x — e~ 1+ ~^) is the exact solution (see Figure S]), with a = IE — 03 in 
order to generate a strong boundary layer. Here, the same refinement procedure than in 
section 15.21 is used. 

Figure [5] shows some of the meshes obtained during the local refinement process. More- 
over, Table [5] displays the corresponding quantitative results. Provided that the boundary 
layer mesh resolution is sufficient, the same behaviours than in the previous tests can also 
be observed : the error goes towards zero as theoretically expected, the effectivity index 
remains almost constant, the superconvergence property of \\a~ l l 2 (Guh — aVu)|| occurs, 
and the mesh is automatically refined around the boundary layer. 
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Figure 5: Mesh levels 1, 5 and 16, non homogeneous boundary-layer case. 
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N 


\\u - u h \\ DG;h 


rv 

^ v error 


V 


Eff 


\\a-^ 2 (Gu h - aVu)\\ 


nv 

^ v recov 


156 


5.90E+00 




9.11E+00 


1.54 


1.46E+00 




766 


1.32E+00 


1.88 


2.52E+00 


1.90 


8.84E-01 


0.63 


3033 


5.76E-01 


1.21 


1.24E+00 


2.16 


3.67E-01 


1.27 


12881 


2.53E-01 


1.13 


5.34E-01 


2.10 


1.19E-01 


1.58 


50496 


1.26E-01 


1.02 


2.66E-01 


2.11 


5.02E-02 


1.26 



Table 5: Boundary-layer case, e = IE — 02, 7 = 250, 7 aj6 = e. 
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